logo

Notes and Conventions

Notes

  • This is an example R Markdown report to demonstrate the workflow for building a “data memo” or report object in parallel with the conduct of a RsNLME analysis.
  • We use html output and “tabsets” within the RMarkdown document to show how RsNLME can be used together with R markdown to produce very dense, interactive data summaries that capture a large amount of data in well organized, easily shared delivery mechanism.

Conventions used in this document

  • We’ve used blue text boxes to indicate descriptive information
  • We’ve used yellow text boxes to indicate results
  • Green text boxes provide key messages from the analysis results
  • Example: This is a key message statement
# Analysis ID     : Certara.R School Lesson 4 - Modeling with RsNLME
# Purpose         : Introduce R Markdown 
# Other Info      : 
library(flextable)
library(plotly)

Objectives

The purpose of this data memo is to provide a summary of the base population pharmacokinetic model fit for study XYZ

Study Design

  • 150 healthy male and female subjects were enrolled in this study
  • Subjects were administered a single 5000, 10000, or 20000 ug dose of Drug XYZ at time=0
  • Blood samples were obtained at 0, 0.25, 0.5, 1, 2, 3, 4, 6, 8, 12, 16 and 24 hours after dosing
  • Covariates in the dataset include baseline subject age, weight, sex, and race

Methods

Overview of the data

Summary

  • Elimination profile appears to be mono-exponential
  • Body weight and/or gender appear to influence PK



Demography summary

dm_table <- readRDS(file='../Lesson_4/dmtable.RDS')
dm_table
Characteristic N = 1501
AGE 50.5 (38.0 - 61.0)
WT 71 (55 - 94)
SEX
F 45 (30%)
M 105 (70%)
RACE
ASIAN 21 (14%)
BLACK 36 (24%)
OTHER 33 (22%)
WHITE 60 (40%)
1 Mean (Range); n (%)



Mean Concentration Summary

mean_conc_table <- readRDS(file='../Lesson_4/mean_conc_table.RDS')
mean_conc_table



PK profiles Linear

pkplot <- readRDS(file='../Lesson_4/pkplot.RDS')

ggplotly(pkplot)



PK profiles Log

  • Disposition appears mono-exponential
pkplotlog <- readRDS(file='../Lesson_4/pkplotlog.RDS')
pkplotlog



Facet Gender and Dose Linear

  • Possible gender difference in PK with observed concentrations higher in females
pkplot_dosesex <- readRDS(file='../Lesson_4/pkplot_dosesex.RDS')
pkplot_dosesex



Facet Race and Dose Linear

  • PK appears similar across race
pkplot_doserace <- readRDS(file='../Lesson_4/pkplot_doserace.RDS')
pkplot_doserace



Facet WT Linear

  • Possible WT influence on PK with larger subjects having lower exposures
pkplot_dosewt <- readRDS(file='../Lesson_4/pkplot_dosewt.RDS')
pkplot_dosewt



Facet Age Linear

  • PK appears similar across age groups. Note limited age range (38 to 61 y).
pkplot_doseage <- readRDS(file='../Lesson_4/pkplot_doseage.RDS')
pkplot_doseage

Base Model Development

Summary

  • A 1 compartment model with first order oral input and elimination was used to fit the data
  • Convergence was achieved with successful estimation of parameter standard errors



Model Description

  • A 1 compartment model with first order oral input and elimination was used to fit the data
  • Model parameters absorption rate constant (Ka), central clearance (Cl), and central volume (V)
  • Between subject variability terms on Ka, Cl and V
  • Additive residual error model

Structural PK Model

## View the updated model 
basemod<-readRDS("../Lesson_4/basemod.RDS")
print(basemod)
## Loading required package: Certara.RsNLME
## 
##  Model Overview 
##  ------------------------------------------- 
## Is population     :  TRUE
## Model Type        :  PK
## 
##  PK 
##  ------------------------------------------- 
## Parameterization  :  Clearance
## Absorption        :  FirstOrder
## Num Compartments  :  1
## Dose Tlag?        :  FALSE
## Elimination Comp ?:  FALSE
## Infusion Allowed ?:  FALSE
## Sequential        :  FALSE
## Freeze PK         :  FALSE
## 
##  PML 
##  ------------------------------------------- 
## test(){
##     cfMicro(A1,Cl/V, first = (Aa = Ka))
##     dosepoint(Aa)
##     C = A1 / V
##     error(CEps=5)
##     observe(CObs=C + CEps)
##     stparm(Ka = tvKa * exp(nKa))
##     stparm(V = tvV * exp(nV))
##     stparm(Cl = (tvCl )  * exp(nCl))
##     fcovariate(AGE)
##     fcovariate(WT)
##     fcovariate(SEX())
##     fcovariate(RACE())
##     fcovariate(DOSEGRP)
##     fixef( tvKa = c(,1.5,))
##     fixef( tvV = c(,80,))
##     fixef( tvCl = c(,9,))
##     ranef(diag(nKa) = c(0.1),diag(nV,nCl) = c(0.1,0.1))
## }
## 
##  Structural Parameters 
##  ------------------------------------------- 
##  Ka V Cl
##  ------------------------------------------- 
## Observations:
## Observation Name :  CObs
## Effect Name      :  C
## Epsilon Name     :  CEps
## Epsilon Type     :  Additive
## Epsilon frozen   :  FALSE
## is BQL           :  FALSE
##  ------------------------------------------- 
##  Column Mappings 
##  ------------------------------------------- 
## Model Variable Name : Data Column name
## id                  : ID
## time                : TIME
## Aa                  : AMT
## AGE                 : AGE
## WT                  : WT
## SEX                 : SEX( M=0 F=1 )
## RACE                : RACE( WHITE=0 BLACK=1 ASIAN=2 OTHER=3 )
## DOSEGRP             : DOSEGRP
## CObs                : CONC

Model Fit Summary

  • Convergence was achieved with successful estimation of parameter standard errors

Additive Error Model

## View the updated model 
basemodfit<-readRDS("../Lesson_4/basemodfit.RDS")
print(basemodfit)
## $Overall
##    Scenario RetCode    LogLik     -2LL      AIC      BIC nParm nObs nSub
## 1: WorkFlow       1 -7187.915 14375.83 14389.83 14427.69     7 1650  150
##    EpsShrinkage Condition
## 1:       0.1073  73.23812
## 
## $theta
##    #Scenario Parameter  Estimate Units     Stderr      CV%    2.5%CI   97.5%CI
## 1:  WorkFlow      tvKa  1.197627    NA 0.03944046 3.293217  1.120268  1.274986
## 2:  WorkFlow       tvV 80.763973    NA 2.84584105 3.523652 75.182115 86.345831
## 3:  WorkFlow      tvCl  7.754275    NA 0.41594592 5.364085  6.938435  8.570115
## 4:  WorkFlow      CEps 13.793006    NA 0.77512135 5.619670 12.272676 15.313335
##    Var.Inf.factor
## 1:             NA
## 2:             NA
## 3:             NA
## 4:             NA
## 
## $omega
##    Label        nKa        nV       nCl
## 1:   nKa 0.07217351 0.0000000 0.0000000
## 2:    nV 0.00000000 0.1695824 0.0000000
## 3:   nCl 0.00000000 0.0000000 0.3959583
## 
## $omega_Correlation
##    Label nKa nV nCl
## 1:   nKa   1  0   0
## 2:    nV   0  1   0
## 3:   nCl   0  0   1
## 
## $Eta_Shrinkage
##        Label       nKa         nV        nCl
## 1: Shrinkage 0.3223916 0.05845163 0.06174757
## 
## $omega_stderr
##    Label        nKa         nV        nCl
## 1:   nKa 0.01528426 0.00000000 0.00000000
## 2:    nV 0.00000000 0.02313111 0.00000000
## 3:   nCl 0.00000000 0.00000000 0.04615598

Multiplicative Error Model

## View the updated model 
basemodNewfit<-readRDS("../Lesson_6/basemodNewfit.RDS")
print(basemodNewfit)
## $Overall
##    Scenario RetCode    LogLik     -2LL      AIC   BIC nParm nObs nSub
## 1: WorkFlow       1 -6744.072 13488.14 13502.14 13540     7 1650  150
##    EpsShrinkage Condition
## 1:      0.12183  87.52661
## 
## $theta
##    #Scenario Parameter  Estimate Units      Stderr      CV%     2.5%CI
## 1:  WorkFlow      tvKa  1.176146    NA 0.034678819 2.948513  1.1081268
## 2:  WorkFlow       tvV 81.436356    NA 2.971941874 3.649404 75.6071625
## 3:  WorkFlow      tvCl  7.947471    NA 0.418135465 5.261240  7.1273360
## 4:  WorkFlow      CEps  0.206732    NA 0.009885775 4.781929  0.1873419
##      97.5%CI Var.Inf.factor
## 1:  1.244165             NA
## 2: 87.265549             NA
## 3:  8.767605             NA
## 4:  0.226122             NA
## 
## $omega
##    Label        nKa        nV       nCl
## 1:   nKa 0.07821959 0.0000000 0.0000000
## 2:    nV 0.00000000 0.1783194 0.0000000
## 3:   nCl 0.00000000 0.0000000 0.3698158
## 
## $omega_Correlation
##    Label nKa nV nCl
## 1:   nKa   1  0   0
## 2:    nV   0  1   0
## 3:   nCl   0  0   1
## 
## $Eta_Shrinkage
##        Label       nKa        nV        nCl
## 1: Shrinkage 0.2143318 0.0331527 0.02153311
## 
## $omega_stderr
##    Label        nKa         nV        nCl
## 1:   nKa 0.01487545 0.00000000 0.00000000
## 2:    nV 0.00000000 0.02313373 0.00000000
## 3:   nCl 0.00000000 0.00000000 0.03852164

Base Model Diagnostics

Additive Residual Error Model

Structural Model

  • A multiplicative residual error model was selected in place of
  • the base model additive error model
DV vs PRED
basedvpred <- readRDS(file='../Lesson_6/basedvpred.RDS')
basedvpred

DV vs IPRED
basedvipred <- readRDS(file='../Lesson_6/basedvipred.RDS')
basedvipred

CWRES vs TIME
baseresidv <- readRDS(file='../Lesson_6/baseresidv.RDS')
baseresidv

CWRES vs IPRED
baserespred <- readRDS(file='../Lesson_6/baserespred.RDS')
baserespred

absCWRES vs IPRED
basearespred <- readRDS(file='../Lesson_6/basearespred.RDS')
basearespred

Covariate Model

  • CL appears to be influenced by subject bodyweight
Parameter vs WT
baseprmcovwt <- readRDS(file='../Lesson_6/baseprmcovwt.RDS')
baseprmcovwt

ETA vs WT
baseetacovwt <- readRDS(file='../Lesson_6/baseetacovwt.RDS')
baseetacovwt

Multiplicative Residual Error Model

Structural Model

  • A multiplicative residual error model was selected in place of
  • the base model additive error model
DV vs PRED
basenewdvpred <- readRDS(file='../Lesson_6/basenewdvpred.RDS')
basenewdvpred

DV vs IPRED
basenewdvipred <- readRDS(file='../Lesson_6/basenewdvipred.RDS')
basenewdvipred

CWRES vs TIME
basenewresidv <- readRDS(file='../Lesson_6/basenewresidv.RDS')
basenewresidv

CWRES vs IPRED
basenewrespred <- readRDS(file='../Lesson_6/basenewrespred.RDS')
basenewrespred

absCWRES vs IPRED
basenewarespred <- readRDS(file='../Lesson_6/basenewarespred.RDS')
basenewarespred

Covariate Model

  • CL appears to be influenced by subject bodyweight
Parameter vs WT
basenewprmcovwt <- readRDS(file='../Lesson_6/basenewprmcovwt.RDS')
basenewprmcovwt

ETA vs WT
basenewetacovwt <- readRDS(file='../Lesson_6/basenewetacovwt.RDS')
basenewetacovwt

ETAs vs All Covariates
basenewranparcov <- readRDS(file='../Lesson_6/basenewranparcov.RDS')
basenewranparcov
## [[1]]

Covariate Models

Summary

  • Significant covariates in single covariate model runs include…



Run1WTonCL

Model Fit Summary

run1WTonCLfit<-readRDS("../Lesson_6/run1WTonCLfit.RDS")
print(run1WTonCLfit)
## $Overall
##    Scenario RetCode    LogLik     -2LL      AIC      BIC nParm nObs nSub
## 1: WorkFlow       1 -6717.738 13435.48 13451.48 13494.74     8 1650  150
##    EpsShrinkage Condition
## 1:      0.12186  87.53302
## 
## $theta
##    #Scenario Parameter   Estimate Units      Stderr       CV%     2.5%CI
## 1:  WorkFlow      tvKa  1.1775716    NA 0.034644115  2.941997  1.1096203
## 2:  WorkFlow       tvV 81.5451627    NA 2.969473676  3.641508 75.7208080
## 3:  WorkFlow      tvCl  7.7866374    NA 0.346324513  4.447677  7.1073531
## 4:  WorkFlow    dCldWT  2.7622329    NA 0.335075232 12.130593  2.1050131
## 5:  WorkFlow      CEps  0.2066357    NA 0.009884458  4.783518  0.1872483
##       97.5%CI Var.Inf.factor
## 1:  1.2455229             NA
## 2: 87.3695174             NA
## 3:  8.4659217             NA
## 4:  3.4194527             NA
## 5:  0.2260232             NA
## 
## $omega
##    Label        nKa        nV       nCl
## 1:   nKa 0.08048666 0.0000000 0.0000000
## 2:    nV 0.00000000 0.1781333 0.0000000
## 3:   nCl 0.00000000 0.0000000 0.2600164
## 
## $omega_Correlation
##    Label nKa nV nCl
## 1:   nKa   1  0   0
## 2:    nV   0  1   0
## 3:   nCl   0  0   1
## 
## $Eta_Shrinkage
##        Label       nKa         nV        nCl
## 1: Shrinkage 0.2186306 0.03389274 0.02339663
## 
## $omega_stderr
##    Label        nKa         nV        nCl
## 1:   nKa 0.01538703 0.00000000 0.00000000
## 2:    nV 0.00000000 0.02313116 0.00000000
## 3:   nCl 0.00000000 0.00000000 0.02590652

Structural Model

DV vs PRED

run1dvpred <- readRDS(file='../Lesson_6/run1dvpred.RDS')
run1dvpred

DV vs IPRED

run1dvipred <- readRDS(file='../Lesson_6/run1dvipred.RDS')
run1dvipred

CWRES vs TIME

run1residv <- readRDS(file='../Lesson_6/run1residv.RDS')
run1residv

CWRES vs IPRED

run1respred <- readRDS(file='../Lesson_6/run1respred.RDS')
run1respred

Covariate Model

ETAs vs All Covariates

run1ranparcov <- readRDS(file='../Lesson_6/run1ranparcov.RDS')
run1ranparcov
## [[1]]

Run2SEXonCL

Model Fit Summary

run2SEXonCLfit<-readRDS("../Lesson_6/run2SEXonCLfit.RDS")
print(run2SEXonCLfit)
## $Overall
##    Scenario RetCode    LogLik     -2LL      AIC      BIC nParm nObs nSub
## 1: WorkFlow       1 -6742.831 13485.66 13501.66 13544.93     8 1650  150
##    EpsShrinkage Condition
## 1:      0.12228  87.79989
## 
## $theta
##    #Scenario Parameter   Estimate Units      Stderr        CV%     2.5%CI
## 1:  WorkFlow      tvKa  1.1756594    NA 0.034666460   2.948682  1.1076642
## 2:  WorkFlow       tvV 81.4338477    NA 2.947980415   3.620092 75.6516501
## 3:  WorkFlow      tvCl  8.3413248    NA 0.449809304   5.392540  7.4590644
## 4:  WorkFlow  dCldSEX1 -0.1748848    NA 0.059423076 -33.978420 -0.2914378
## 5:  WorkFlow      CEps  0.2066139    NA 0.009851288   4.767969  0.1872915
##        97.5%CI Var.Inf.factor
## 1:  1.24365450             NA
## 2: 87.21604528             NA
## 3:  9.22358514             NA
## 4: -0.05833175             NA
## 5:  0.22593635             NA
## 
## $omega
##    Label        nKa        nV       nCl
## 1:   nKa 0.08148876 0.0000000 0.0000000
## 2:    nV 0.00000000 0.1779236 0.0000000
## 3:   nCl 0.00000000 0.0000000 0.3632091
## 
## $omega_Correlation
##    Label nKa nV nCl
## 1:   nKa   1  0   0
## 2:    nV   0  1   0
## 3:   nCl   0  0   1
## 
## $Eta_Shrinkage
##        Label       nKa         nV        nCl
## 1: Shrinkage 0.2179584 0.03227745 0.02093563
## 
## $omega_stderr
##    Label        nKa         nV        nCl
## 1:   nKa 0.01495687 0.00000000 0.00000000
## 2:    nV 0.00000000 0.02277476 0.00000000
## 3:   nCl 0.00000000 0.00000000 0.03709178

Structural Model

DV vs PRED

run2dvpred <- readRDS(file='../Lesson_6/run2dvpred.RDS')
run2dvpred

DV vs IPRED

run2dvipred <- readRDS(file='../Lesson_6/run2dvipred.RDS')
run2dvipred

CWRES vs TIME

run2residv <- readRDS(file='../Lesson_6/run2residv.RDS')
run2residv

CWRES vs IPRED

run2respred <- readRDS(file='../Lesson_6/run2respred.RDS')
run2respred

Covariate Model

ETAs vs All Covariates

run2ranparcov <- readRDS(file='../Lesson_6/run2ranparcov.RDS')
run2ranparcov
## [[1]]

Run3WTonV

Model Fit Summary

run3WTonVfit<-readRDS("../Lesson_6/run3WTonVfit.RDS")
print(run3WTonVfit)
## $Overall
##    Scenario RetCode    LogLik  -2LL   AIC      BIC nParm nObs nSub EpsShrinkage
## 1: WorkFlow       1 -6744.001 13488 13504 13547.27     8 1650  150      0.12163
##    Condition
## 1:  86.58998
## 
## $theta
##    #Scenario Parameter    Estimate Units      Stderr       CV%     2.5%CI
## 1:  WorkFlow      tvKa  1.17641526    NA 0.034686689  2.948507  1.1083805
## 2:  WorkFlow       tvV 81.39156156    NA 2.937142300  3.608657 75.6306219
## 3:  WorkFlow      tvCl  7.93844189    NA 0.416285225  5.243916  7.1219360
## 4:  WorkFlow     dVdWT  0.09994762    NA 0.094422862 94.472342 -0.0852543
## 5:  WorkFlow      CEps  0.20651569    NA 0.009874919  4.781680  0.1871469
##       97.5%CI Var.Inf.factor
## 1:  1.2444501             NA
## 2: 87.1525012             NA
## 3:  8.7549478             NA
## 4:  0.2851495             NA
## 5:  0.2258844             NA
## 
## $omega
##    Label        nKa        nV       nCl
## 1:   nKa 0.08045648 0.0000000 0.0000000
## 2:    nV 0.00000000 0.1787547 0.0000000
## 3:   nCl 0.00000000 0.0000000 0.3701807
## 
## $omega_Correlation
##    Label nKa nV nCl
## 1:   nKa   1  0   0
## 2:    nV   0  1   0
## 3:   nCl   0  0   1
## 
## $Eta_Shrinkage
##        Label       nKa         nV        nCl
## 1: Shrinkage 0.2166663 0.03497188 0.02198844
## 
## $omega_stderr
##    Label        nKa         nV        nCl
## 1:   nKa 0.01552365 0.00000000 0.00000000
## 2:    nV 0.00000000 0.02328333 0.00000000
## 3:   nCl 0.00000000 0.00000000 0.03884055

Structural Model

DV vs PRED

run3dvpred <- readRDS(file='../Lesson_6/run3dvpred.RDS')
run3dvpred

DV vs IPRED

run3dvipred <- readRDS(file='../Lesson_6/run3dvipred.RDS')
run3dvipred

CWRES vs TIME

run3residv <- readRDS(file='../Lesson_6/run3residv.RDS')
run3residv

CWRES vs IPRED

run3respred <- readRDS(file='../Lesson_6/run3respred.RDS')
run3respred

Covariate Model

ETAs vs All Covariates

run3ranparcov <- readRDS(file='../Lesson_6/run3ranparcov.RDS')
run3ranparcov
## [[1]]

Simulation

Summary

  • The final model was used to simulate a multiple dose regime of 5000, 10000, or 20000 ug q6h * 2 days …



Simulation Output

finalsimq6hbywtgrp<-readRDS("../Lesson_7/finalsimq6hbywtgrp.RDS")
print(finalsimq6hbywtgrp)

Session Info

sessionInfo()
R version 4.2.1 (2022-06-23 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19044)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.utf8 
[2] LC_CTYPE=English_United States.utf8   
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C                          
[5] LC_TIME=English_United States.utf8    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] Certara.RsNLME_1.1.0 plotly_4.10.0        ggplot2_3.3.6       
[4] flextable_0.7.3     

loaded via a namespace (and not attached):
 [1] nlme_3.1-157             RColorBrewer_1.1-3       httr_1.4.4              
 [4] tools_4.2.1              backports_1.4.1          bslib_0.4.0             
 [7] utf8_1.2.2               R6_2.5.1                 DBI_1.1.3               
[10] lazyeval_0.2.2           mgcv_1.8-40              colorspace_2.0-3        
[13] withr_2.5.0              gridExtra_2.3            GGally_2.1.2            
[16] tidyselect_1.1.2         curl_4.3.2               compiler_4.2.1          
[19] Certara.Xpose.NLME_1.1.0 cli_3.3.0                gt_0.6.0                
[22] xml2_1.3.3               shinyjs_2.1.0            officer_0.4.3           
[25] shinymaterial_1.2.0      labeling_0.4.2           sass_0.4.2              
[28] scales_1.2.1             checkmate_2.1.0          commonmark_1.8.0        
[31] systemfonts_1.0.4        stringr_1.4.1            digest_0.6.29           
[34] rmarkdown_2.16           katex_1.4.0              base64enc_0.1-3         
[37] Certara.NLME8_1.2.1      pkgconfig_2.0.3          htmltools_0.5.3         
[40] fastmap_1.1.0            highr_0.9                htmlwidgets_1.5.4       
[43] rlang_1.0.4              rstudioapi_0.13          shiny_1.7.2             
[46] jquerylib_0.1.4          farver_2.1.1             generics_0.1.3          
[49] jsonlite_1.8.0           crosstalk_1.2.0          dplyr_1.0.9             
[52] zip_2.2.0                magrittr_2.0.3           Matrix_1.4-1            
[55] Rcpp_1.0.9               munsell_0.5.0            fansi_1.0.3             
[58] gdtools_0.2.4            lifecycle_1.0.1          stringi_1.7.8           
[61] yaml_2.3.5               gtsummary_1.6.1          MASS_7.3-57             
[64] plyr_1.8.7               grid_4.2.1               promises_1.2.0.1        
[67] crayon_1.5.1             egg_0.4.5                lattice_0.20-45         
[70] splines_4.2.1            knitr_1.40               pillar_1.8.1            
[73] uuid_1.1-0               xslt_1.4.3               glue_1.6.2              
[76] evaluate_0.16            V8_4.2.1                 data.table_1.14.2       
[79] broom.helpers_1.8.0      vctrs_0.4.1              tweenr_1.0.2            
[82] httpuv_1.6.5             gtable_0.3.0             purrr_0.3.4             
[85] polyclip_1.10-0          tidyr_1.2.0              reshape_0.8.9           
[88] xpose_0.4.13             assertthat_0.2.1         cachem_1.0.6            
[91] xfun_0.32                ggforce_0.3.3            mime_0.12               
[94] xtable_1.8-4             later_1.3.0              equatags_0.2.0          
[97] viridisLite_0.4.1        tibble_3.1.8             ellipsis_0.3.2